clear all; close all; clear classes; clc;

%% Palik's Ag, p.353
wvlen1 = 1e-1 * [1.240 1.550 2.066 2.480 2.755 3.099 3.444 3.542 4.133 4.959 6.199 8.265 12.40 15.50 16.75 20.66 24.80 30.99 33.51 35.42 41.33 44.28 47.68 51.66 56.35 61.99 68.88 72.93 77.49 82.65 88.56 95.37 103.3 112.7 118.1 124.0 130.5 137.8 145.9 155.0 158.9 163.1 167.5 172.2 177.1 182.3 187.8 193.7 200.0 206.6 213.8 221.4 229.6 233.9];
n1 = [0.9990 1.000 1.002 1.001 1.001 1.003 1.003 1.003 1.001 1.001 1.001 0.999 0.998 0.997 0.998 0.997 0.997 0.998 1.007 1.001 0.996 0.994 0.992 0.989 0.987 0.982 0.978 0.973 0.968 0.963 0.955 0.943 0.926 0.902 0.890 0.876 0.861 0.849 0.846 0.848 0.851 0.853 0.853 0.851 0.859 0.871 0.883 0.885 0.881 0.873 0.868 0.871 0.885 0.890];

k11 = 1e-6 * [1.29 2.92 8.36];
k12 = 1e-5 * [1.52 2.21 3.57 3.24 1.72 1.43 2.40 4.72];
k13 = 1e-4 * [1.21 5.23];
k14 = 1e-3 * [1.08 1.42 2.69 4.57 8.82 8.23 1.68 2.32 2.70 3.11 3.50 3.45 3.77 4.31 4.42 4.38 4.07 3.54 3.66 6.22];
k15 = 1e-2 * [1.72 2.61 3.82 5.65 8.24];
k16 = [0.111 0.139 0.149 0.158 0.166 0.181 0.197 0.208 0.211 0.210 0.211 0.221 0.237 0.259 0.274 0.275];

k1 = [k11 k12 k13 k14 k15 k16];


%% Palik's Ag, p.354
wvlen_n2 = [
    238.4	0.889 
    243.1	0.886 
    248.0	0.884 
    253.0	0.885 
    258.3	0.888 
    263.8	0.893 
    269.5	0.896 
    281.8	0.899 
    295.2	0.897 
    309.9	0.896 
    317.9	0.895 
    326.3	0.892 
    335.1	0.885 
    344.4	0.876 
    354.2	0.865 
    364.6	0.879 
    375.7	0.899 
    387.4	0.921 
    399.9	0.932 
    413.2	0.931 
    427.5	0.919 
    435.0	0.911 
    442.8	0.902 
    450.8	0.851 
    459.2	0.886
    467.9	0.886 
    476.9	0.887 
    486.2	0.891 
    506.1	0.906 
    527.6	0.943 
    539.1	0.980 
    551.0	1.029 
    563.6	1.075 
    576.7	1.104 
    590.4	1.109 
    604.8	1.105 
    619.9	1.098 
    652.6	1.096 
    688.8	1.101 
    729.3	1.121 
    774.9	1.167 
    799.9	1.202 
    826.6	1.243 
    855.1	1.281 
    885.6	1.308 
    918.4	1.322 
    953.7	1.322 
    991.9	1.315
    ];

wvlen2 = 1e-1 * wvlen_n2(:,1).';
n2 = wvlen_n2(:,2).';

k2 = [0.276 0.281 0.290 0.300 0.309 0.317 0.323 0.334 0.349 0.368 0.378 0.388 0.400 0.418 0.454 0.489 0.514 0.528 0.534 0.541 0.557 0.572 0.590 0.616 0.650 0.421 0.444 0.468 0.522 0.584 0.614 0.624 0.615 0.584 0.561 0.549 0.548 0.565 0.593 0.635 0.676 0.691 0.693 0.679 0.656 0.628 0.603 0.588];


%% Palik's Ag, p.354, micron section
wvlen_n3 = [
    0.1025	1.308
    0.1051	1.304
    ];

wvlen3 = 1e3 * wvlen_n3(:,1).';
n3 = wvlen_n3(:,2).';

k3 = [0.581 0.578];


%% Palik's Ag, p.355
wvlen4 = 1e3 * [0.1078 0.1107 0.1137 0.1170 0.1204 0.1240 0.1278 0.1348 0.1378 0.1409 0.1442 0.1476 0.1512 0.1550 0.1590 0.1631 0.1675 0.1722 0.1771 0.1823 0.1879 0.1937 0.1968 0.2000 0.2033 0.2066 0.2138 0.2214 0.2296 0.2384 0.2480 0.2530 0.2583 0.2638 0.2695 0.2755 0.2818 0.2883 0.2952 0.2988 0.3024 0.3061 0.3100 0.3115 0.3139 0.3155 0.3179 0.3195 0.3220 0.3237 0.3263 0.3306 0.3324 0.3351];
n4 = [1.300 1.293 1.280 1.265 1.252 1.241 1.229 1.182 1.149 1.112 1.073 1.032 0.993 0.962 0.940 0.935 0.936 0.942 0.953 0.969 0.995 1.028 1.048 1.072 1.098 1.125 1.173 1.208 1.238 1.265 1.298 1.320 1.343 1.372 1.404 1.441 1.476 1.502 1.519 1.522 1.496 1.432 1.323 1.246 1.149 1.044 0.932 0.815 0.708 0.616 0.526 0.371 0.321 0.294];
k4 = [0.573 0.566 0.560 0.560 0.564 0.568 0.566 0.550 0.552 0.563 0.581 0.610 0.653 0.706 0.770 0.832 0.892 0.951 1.01 1.07 1.13 1.18 1.21 1.24 1.26 1.27 1.29 1.30 1.31 1.33 1.35 1.35 1.35 1.35 1.33 1.31 1.26 1.19 1.08 0.992 0.882 0.766 0.647 0.586 0.540 0.514 0.504 0.526 0.565 0.609 0.663 0.813 0.902 0.986];


%% Palik's Ag, p.356
wvlen5 = 1e3 * [0.3397 0.3444 0.3542 0.3647 0.3757 0.3875 0.4000 0.4133 0.4275 0.4428 0.4592 0.4769 0.4959 0.5166 0.5391 0.5636 0.5904 0.6199 0.6526 0.6888 0.7293 0.7749 0.8266 0.8856 0.9537 1.033 1.127 1.240 1.265 1.291 1.305 1.319 1.348 1.378 1.409 1.442 1.459 1.476 1.512 1.550 1.590 1.631 1.653 1.675 1.722 1.771 1.823 1.879 1.907 1.937 2.000];
n5 = [0.259 0.238 0.209 0.186 0.200 0.192 0.173 0.173 0.160 0.157 0.144 0.132 0.130 0.130 0.129 0.120 0.121 0.131 0.140 0.140 0.148 0.143 0.145 0.163 0.198 0.226 0.251 0.329 0.375 0.383 0.358 0.392 0.401 0.411 0.421 0.431 0.446 0.442 0.455 0.469 0.485 0.501 0.624 0.519 0.537 0.557 0.578 0.600 0.873 0.624 0.650];
k5 = [1.12 1.24 1.44 1.61 1.67 1.81 1.95 2.11 2.26 2.40 2.56 2.72 2.88 3.07 3.25 3.45 3.66 3.88 4.15 4.44 4.74 5.09 5.50 5.95 6.43 6.99 7.67 8.49 7.78 7.92 8.95 8.06 8.21 8.37 8.37 8.70 10.1 8.88 9.08 9.32 9.57 9.84 11.5 10.1 10.4 10.7 11.1 11.4 13.3 11.8 12.2];


%% Palik's Ag, p.357
wvlen6 = 1e3 * [2.066 2.138 2.214 2.296 2.384 2.480 2.583 2.695 2.818 2.952 3.100 3.263 3.444 3.647 3.875 4.133 4.428 4.769 5.166 5.636 6.199 6.526 6.888 7.293 7.749 8.266 8.856 9.537 9.919];
n6 = [0.668 0.729 0.774 0.823 0.878 0.939 1.007 1.083 1.168 1.265 1.387 1.536 1.710 1.915 2.160 2.446 2.786 3.202 3.732 4.425 5.355 5.960 6.670 7.461 8.376 9.441 10.69 12.21 13.11];
k6 = [12.6 13.0 13.5 14.0 14.5 15.1 15.7 16.4 17.1 17.9 18.8 19.8 20.9 22.1 23.5 25.1 26.9 29.0 31.3 34.0 37.0 38.6 40.4 42.5 44.8 47.1 49.4 52.2 53.7];

%% Collect all.
wvlen = [wvlen1 wvlen2 wvlen3 wvlen4 wvlen5 wvlen6].';
n = [n1 n2 n3 n4 n5 n6].';
k = [k1 k2 k3 k4 k5 k6].';

%% Convert the wavelengths to the photon energies.
eV = PhysC.h * PhysC.c0 * 1e9 ./ wvlen;

%% Reverse the data order.
eV = eV(end:-1:1);
n = n(end:-1:1);
k = k(end:-1:1);
wvlen = wvlen(end:-1:1);

%% Calculate the permittivity from n and k following the exp(+i w t) time dependence.
eps = (n - 1i*k).^2;

%% Plot parameters.  Compare with Fig.9 on p.352 of Palik.
nk_wvlen = 1;
eps_eV = 2;
eps_wvlen = 3;
plotstyle = nk_wvlen;
switch plotstyle
    case nk_wvlen  % plot n and k
        loglog(wvlen, n, 'o-', wvlen, k, 'o-')
        %plot(wvlen, n, wvlen, k)
        legend('n', 'k', 'Location', 'SouthEast');
        xlabel 'wavelength (nm)'
        %axis([1e2 1e4 1e-2 1e2])
    case eps_eV  % plot real(eps) and -imag(eps)
        plot(eV, real(eps), 'o-', eV, -imag(eps), 'o-')
        legend('\epsilon_1', '\epsilon_2', 'Location', 'SouthEast');
        xlabel 'Photon Energy (eV)'
        axis([0.5 6.5 -7 7]);
    case eps_wvlen
        plot(wvlen, real(eps), 'o-', wvlen, -imag(eps), 'o-')
        legend('\epsilon_1', '\epsilon_2', 'Location', 'SouthEast');
        xlabel 'wavelength (nm)'
        %axis([0 1e2 -1e1 1e1])
end
title(mfilename)

%% Save data.
save(mfilename, 'eV', 'n', 'k');
